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Abstract. We consider the task of intrinsic control system identification for quantum devices. 

The problem of experimental determination of subspaco confinement is considered, and simple 
general strategies for full Hamiltonian identification and decoherence characterization of a 
controlled two-level system are presented. 



1. Introduction 

Advances in nano-fabrication are increasingly enabling us to create nano-scale devices that 
exhibit non-classical or quantum-mechanical behaviour. Such quantum devices are of great 
interest as they may pave the way for a new generation of quantum technology with various 
applications from quantum metrology to quantum information processing. However, to create 
quantum devices that perform useful functions, we must be able to understand their behaviour, 
and have effective means to controllably manipulate it. Analysis of system dynamics and the 
design of effective control strategies is almost impossible without the availability of sufficiently 
accurate mathematical models of the device. While these models should capture the essential 
features of the device, to be useful, they must also be computationally tractable, and preferably 
as simple as possible. 

There are different approaches to deal with this problem. One, which we shall refer to 

as a first-principles approach, involves constructing model Hamiltonians based on reasonable 
assumptions about the relevant physical processes governing the behaviour of the system, 
making various simplifications, and solving the Schrodingcr equation in some form, usually using 
numerical techniques such as finite element or functional expansion methods. The empirical 
approach, on the other hand, starts with experimental data and observations to construct a 
model of the system. In practice both approaches are needed to deal with complex systems. 
First-principle models are crucial to elucidate the fundamental physics that governs a system 
or device, but experimental data is crucial to account for the many unknowns that result in 



Figure 1: Simple charge 'qubit' device: Trench-isolated Silicon double quantum dot molecule 
fabricated using e-beam lithography, controlled by several DC-gate electrodes, with a single 
electron transistor for charge readout (Hitachi Cambridge Labs) (jlj). 

variability of the characteristics of the system, which is often difficult to predict or explain 
theoretically. For instance, 'manufactured' systems such as artificial quantum dot 'atoms' and 
'molecules' vary considerably in size, geometry, internal structure, etc, and sometimes even small 
variations can significantly alter their behaviour. 

In this paper we will focus on constructing mathematical models solely from experimental 
data without any prior knowledge of the specifics of the system, an approach one might call 
black-box system identification. We will restrict ourselves to quantum systems whose essential 
features can be captured by low-dimensional Hilbert space models. For instance, a quantum 
dot molecule such as the Silicon double quantum dot system in Fig. [T] may consist of millions 
of individual atoms with many degrees of freedom, but the observable dynamics may be well 
described by overall charge distribution states in the quantum dot molecule. Taking these charge 
distribution states as basis states for a Hilbert space the observable dynamics of the system 
is governed by a Hamiltonian operator acting on 7i, and possibly additional non-Hermitian 
operators to account for dissipative effects if the system is not completely isolated from its 
environment. If the dimension of the Hilbert space 7i is huge then complete characterization 
of the Hamiltonian and dissipation operators may be a hopeless task, but in certain cases the 
dynamics of interest takes place in low- dimensional Hilbert space or we may in fact wish to 
design a device whose dynamics can be described by a low-dimensional Hilbert space model, as 
is the case in quantum information processing ([1]), where we desire to create systems that act 
as quantum bits, for example. 

The paper is organized as follows. We start with a general description of the basic assumptions 
we make about the system to be modelled and the type of measurements and experiments upon 
which our characterization protocols are based. In section [3] we discuss the problem of subspace 
confinement, i.e., how to experimentally characterize how well the dynamics of system is confined 
to a low-dimensional subspace such as a two-dimensional (qubit) Hilbert space. In section |4] we 
present general protocols for identifying Hamiltonian and decoherence parameters for controlled 
qubit-like systems. We conclude with a discussion of generalizations and open problems. 

2. General Modelling Assumptions 

We start with a generic model assuming only that the state of the system can be represented 
by a density operator (i.e., a positive operator of trace one) p acting on a Hilbert space TL of 
dimension A^, and that its evolution is governed by the quantum Liouville equation (in units 
such that fi = 1): 

|p(t) = -i[i7,p(t)]+L^[p(t)], (1) 

where H is an effective Hamiltonian and L/j a super-operator that accounts for dissipative effects 
due to environmental influences, etc. Initially, all that is known about H and L^) is that H is 
a Hermitian operator on "H, and L^j a completely positive super-operator acting on density 
matrices p{t), although we may make some additional assumptions about the structure of H 
or L/). For instance, we shall generally assumes L/j to be of Kossakowski-Sudarshan-Lindblad 



form daSD: 

k=l 

where the super-operators VlVk] are defined by 

V[A]B = ABA^ - {A^AB + BA^A) /2, (3) 

Vk are (generally non-Hermitian) operators on Ti and are positive real numbers. 

If the system is subject to dynamic control, then the Hamiltonian (and sometimes the 
relaxation operators Vk) are not fixed but dependent on external control fields, which we shall 
denote by f(t). In the simplest case, we may assume a linear dependence of the Hamiltonian on 
the controls 

M 

H[m] =Ho+Yl fm{t)Hm, (4) 
m=l 

where the Hm are fixed Hermitian operators on TC and the reservoir operators Vk are constant. 

The ultimate objective of characterization is to identify the operators if[f(t)] (or Hm) and 
by performing suitable experiments on the system. This task would be simplified if we could 
assume that the system can be prepared in an arbitrary initial state po, and that we can perform 
generalized measurements or projective measurements in arbitrary bases at any time — but these 
requirements are generally unrealistic. For example, our ability to perform measurements on the 
system prior to characterization is limited by the direct readout processes available, and often we 
only have a single, basic sensor such as a single electron transistor (SET) (jS)) providing relatively 
limited information about the charge distribution in a quantum dot molecule, or detection may 
be accomplished via a readout transition that involves coupling one state of the system to a 
fluorescent read-out state, e.g., via a laser, etc. Furthermore, preparation of non-trivial states 
generally depends on knowledge of the operators lf[f(i)] and L^), the very information about the 
system we are trying to obtain. These practical restrictions rule out conventional quantum state 
or process tomography techniques, which presume the ability to measure the system in different 
measurement bases and prepare it in different initial states to obtain sufficient information to 
reconstruct the quantum state or process ([61 [71 (HJ [9l [TO\i . 

In this paper we consider a rather typical experimental scenario, where we are limited to 
measurements of a fixed observable and evolution under a Hamiltonian that can be modified by 
varying certain control settings. For the most part, we restrict ourselves further to piecewise 
constant controls. The only assumption on the measurement process we make is that it can 
be formally represented by some Hermitian operator A with eigenvalues corresponding to 
measurement outcomes A„, i.e., that it has a spectral decomposition of the form 

N 

A = ^Xn\n){n\. (5) 

n=l 

This measurement also serves as initialization of the system as outcome means that the system 
will be left in an eigenstate |n) associated with the eigenvalue A^. If A has N unique eigenvalues, 
i.e., all eigenvalues occur with multiplicity one, then the measurement is sufficient to initialize 
the system in a unique state; if it has degenerate eigenvalues associated with eigenspaces of 
dimension greater > 1 then some measurement outcomes will not determine the state uniquely. 

Following the idea of intrinsic characterization, our objective is to extract information 
about the system without recourse to any external resources, i.e., using no information from 
measurements other than the information provided by the sensors built into the device, and no 



external control fields except the ability to change the settings of the built-in actuators (such 
as variable gate voltages) subject to constraints. Although this requirement of characterization 
relying only on the built-in sensors and actuators may seem excessively restrictive, excluding 
many forms of spectroscopy, for example, it has the advantage of simplicity (no external 
resources required). Moreover, external sensors and actuators may disturb the system, and 
thus characterization of the system in their presence may not yield an accurate picture of the 
dynamics in the absence of the additional apparatus. 

We restrict ourselves here to systems sufficiently weakly coupled to a sufficiently large 
reservoir, whose dynamics can be described by adding a dissipation super-operator of the form ^ 
to the Hamiltonian dynamics of the subsystem of interest. Although many systems can be 
modelled this way, it should be noted that this approach has limitations. For example, the 
dynamics of a subsystem Tis that is strongly is coupled to a finite reservoir TCr, such as single 
spin coupled to several nearby spins, can be very complicated and non-Markovian. Although 
non-Markovian dynamics can in principle be dealt with by allowing time-dependent relaxation 
operators Vk, in such cases it is not always possible to describe the dynamics of the subsystem 
of interest in terms of the Hamiltonian dynamics on the subspace and a set of simple relaxation 
or decoherence operators. Rather, it may become necessary to consider the system Tls + Ti-R 
and characterize its dynamics, described by a Schrodinger equation with a Hamiltonian Hs+r, 
instead to obtain an accurate picture of the subsystem dynamics. The basic ideas of intrinsic 
characterization can be applied to this larger system, and the protocols we will describe can in 
principle be extended to such higher-dimensional systems, although full characterization of the 
system plus reservoir Hamiltonian Hs+r using only the built-in sensors and actuators may not 
be possible. The degree of characterization possible will depend on the size of the Hilbert space 
and the capabilities of the built-in sensors and actuators, e.g., to discriminate and manipulate 
different states of the larger system. Before attempting to identify the system Hamiltonian (and 
decoherence operators), an important first task is therefore estimating the dimension of the 
Hilbert space in which the dynamics takes place. 

3. Characterization of subspace confinement 

A fundamental prerequisite for constructing a Hilbert space model is knowledge of the underlying 
Hilbert space. This is a nontrivial problem as most systems have many degrees of freedom, and 
thus a potentially huge Hilbert space, but effective characterization of the system often depends 
on finding a low dimensional Hilbert space model that captures the essential features of the 
system. Furthermore, in applications such as quantum information processing the elementary 
building blocks are required to have a certain Hilbert space dimension. For example, for a 
quantum-dot molecule to qualify as a qubit, we must be able to isolate a two-dimensional 
subspace of the total Hilbert space, and be able to coherently manipulate states within this 
2D subspace without coupling to states outside the subspace (leakage). This requires several 
characterization steps: 

(i) Isolation of a 2D subspace and characterization of subspace confinement; 

(ii) Characterization of the Hamiltonian dynamics including effect of the actuators; and 

(iii) Characterization of (non-controllable) environmental effects (dissipation). 

The choice of a suitable subspace depends also on the measurement device as the measurement 
must be able to distinguish the basis states. Thus for a potential charge qubit device, for example, 
only a subspace spanned by charge states that can be reliably distinguished by the SET is a 
suitable candidate, and ideally there should be two (or in general N) orthogonal states that 
can be perfectly distinguished by the measurement, so that we have a well-defined reference 
frame (basis) {|0), . . . , |A^)}, and the measurement can be represented by an observable of 
the form ^ with N distinct eigenvalues A„. Furthermore, we would like the readout process 




Figure 2: Leakage characterization assuming three-outcome measurement for system initialized 
in state |0). The continuous curve is the calculated value of pieak(i) = 1 — Tr[n/5(t)], where 
n = |0)(0| + |1)(1| is the projector onto the target subspace, for a 10-level test Hamiltonian. 
The red dots are the estimates of pieak(^) ^-t the discrete times = {k — l)At with = 0.01, 
obtained by averaging over N'^ = 100 simulated experiments for each t}^. As the figure shows, 
the effect of the projection noise makes it difficult to obtain accurate estimates of p\eak{t) for 
any given t = unless N'^ is extremely large. However, the average leakage rate over the time 
interval [0, 100], pieak = 0.11% for both the calculated and measured values. 



to act as a projective measurement, so that each measurement projects the state of the system 
onto one of the measurement basis states The ability to perform projective measurements 
is not a trivial requirement, especially for solid state systems, and is not an absolute necessity as 
characterization protocols can be adapted to weak measurements, but projective measurements 
can in principle be achieved with various sensors such as rf-SETs, and we shall assume here that 
we are working in this regime. 

After a possible subspace has been identified, it is crucial to check that the subspace is 
sufficiently isolated, i.e., that we can reliably initialize the system in a state in this subspace and 
that it remains in this subspace under both free and controlled evolution. This characterization 
of subspace confinement is very important. If the device is a candidate for a qubit, for example, 
it is essential that subspace leakage be much less than the rate of bit or phase flip errors due to 
imperfect control or decoherence as leakage or loss error correction protocols are considerably 
more demanding in terms of complexity and resources required. 

Characterization of subspace confinement depends on the characteristics of the sensors, i.e., 
the type of measurements we can perform. For our charge qubit with SET readout example, if 
the SET can be calibrated to be sufficiently sensitive to enable detection of states outside the 
chosen subspace in addition to being able to discriminate the subspace basis states i.e., if the 
true measurement has (at least) N mutually exclusive outcomes A„, n = 1, . . . , and Aat+i 
if the system state is outside the subspace, then the characterization of subspace confinement is 
relatively easy. If we perform Ag experiments of the form 

(i) Initialize: Measure and record outcome Aq 

(ii) Evolve: Let the system evolve for time t under some fixed Hamiltonian Hf 

(iii) Measure: Repeat measurement and record outcome A;, 

and denote by Na^b the number of times the measurement results for the initial and final 
measurements were Aq and A^, respectively, then, for Ag sufficiently large, the leakage out of the 
subspace is approximately equal to the fraction of experiments for which the first measurement 
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Figure 3: po{t) obtained by averaging over Ne = 100 simulated experiments with a 10-level test 
Hamiltonian and corresponding Fourier spectrum. From the 0th and 1st order peak heights we 
obtain hQ + 2hi = 0.9970, which yields upper and lower bounds for e of « 0.0015, although these 
bounds are statistically uncertain due to projection noise and discretization errors. Nonetheless, 
it can be verified that the actual subspace confinement for the chosen test Hamiltonian is 9' 
and thus the upper bound provides a reasonable estimate for the leakage. 



was A„ with n < N and the second measurement was + 1, i.e. 

Ne-En=lNN+l,n 



Pleak^— • VOj 



Repeating the experiments for different evolution times t and averaging then gives an indication 
of the rate of subspace leakage, e.g., if we estimate pieakC^fc) for = kAt, k = 0, . . . , K, then 

1 / i^max 1 

Weak = / Pleak(t) dt ^ Y^Pleakfa). (7) 

tmax JO A + i 

If the SET (or other measurement process) is not sufficiently sensitive to reliably distinguish 
at least N subspace basis states |n) as well as states outside the subspace, then subspace 
characterization is more challenging. For instance, suppose we have a 2D subspace and a 
sensor that can reliably detect only one state, say |0), a common scenario for many systems 
where readout transitions are used that detect only a single state. In this case the measurement 
outcomes are Aq or ^Xq. If the measurement is projective and the dynamics confined to a 2D 
subspace, we can identify -iAq with outcome Ai but this identification will lead to errors if 
the dynamics is not really confined to a 2D subspace. However, even in this case we can still 
estimate the level of confinement to a 2D subspace under the evolution of a Hamiltonian Hf 
from observable coherent oscillations. For example, let 

Po(t) = |(0|e-^^*|0)p (8) 



be the probability of obtaining outcome Aq when measuring the time-evolved state \^{t)) = 
exp(— itiff )|0). If the system is initialized in the state |0) and the dynamics of the system under 



Hf is perfectly confined to a two- level subspace then conservation of probability implies that the 
heights ho and hi of the 0th and 1st order terms in the Fourier spectrum of po(0 niust satisfy 
hf) + 2hi = 1. We can use the deviation from this equality to bound the subspace leakage e (? ): 



As the upper and lower bounds depend only on the 0th and 1st order Fourier peaks, they 
can usually be easily determined from experimental data although finite resolution due to 
discretization and projection noise will reduce the accuracy of the estimates (? ). Moreover, in 
this case 1 — e only indicates the confinement of the system to some two-level subspace, and we 
must be careful as for different controls f the dynamics under Hf may be confined to different 
two-level subspaces, in which case the system cannot be considered a proper qubit. For instance, 
for the 10-level test Hamiltonian 
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used for the simulated experiments in Figs |2] and [3j the first procedure measures the projection 
onto the subspace Si spanned by the measurement basis states |0) = (1,0, 0, 0, 0, 0, 0, 0, 0, 0, 0)^ 
and |1) = (0,1,0,0,0,0,0,0,0,0,0)"^, while the second procedure measures the confinement of 
the dynamics to the best-fitting 2D subspace 5*2, which for the given Hamiltonian is spanned by 

VI = (1, -0.0001, 0.0004, 0.0034, -0.0012, -0.0002, -0.0005, -0.0004, 0.0003, -0.0005)^, 
V2 = (0, 1, 0.0174, 0.0238, -0.0131, -0.0051, -0.0032, -0.0020, -0.0022, -0.0013)^, 

and comparison of the results confirms that the average confinement 1 — pieak for the subspace 
Si as obtained by the first procedure is less (99.89%) than the confinement for subspace 5*2 
(99.94%). 

4. General characterization protocols for qubit systems 

Once a suitable subspace has been chosen, we can proceed to the second stage of the 
characterization, identification of the Hamiltonian and decoherence operators. The simplest 
type of system we can consider here is a qubit with Hilbert space dimension 2. In this case, any 
non-trivial projective measurement, i.e., any measurement with two distinguishable outcomes 
Ao and Ai can be represented by an observable A = Ao|0)(0| -|- Ai|l)(l|. The eigenstates |0) and 
|1) of ^4 define a basis for the Hilbert space, and we can define the Pauli operators ao, cfx, cry 
and dz with respect to this basis 
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or in the usual matrix notation 
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Figure 4: Bloch sphere representation of dynamics, (a) Arrangement of rotation axes with 
respect to reference rotation axis d^. (b) The evolution of the vector s(t) with sq = (0,0, 1)"^ 
about the axis d, and the projection of s(t) onto the z-axis. 



Taking (without loss of generality) the measurement outcomes to be Aq = 1 and Ai = —1, we 
have ^ = (T^ in this basis. Furthermore, we can expand any Hamiltonian of the system as 



H = ^(doo-Q + d^cTx + dy&y + d;,az) 

— ~2 ^ 4'(^x + sin 6 sin (j)ay + cos 6* a^) . 



(12a) 
(12b) 



The do term can generally be ignored as the identity ctq commutes with all other Pauli matrices 
and exp(— itdoCTo) corresponds to multiplication by a global phase factor, i.e., trivial dynamics. 
Hence, it suffices to determine the real vector d = {dx,dy,dz), or in polar form, the angles 
{uj,9,(j)) to determine the Hamiltonian. For a single Hamiltonian, we can furthermore choose 
the coordinate system such that (p = 0. However, since H depends on control inputs f, the 
parameters co, 9 and (p also depend on the controls f , and we usually need to determine Hf for 
many different f. In this case we can choose a reference Hamiltonian, e.g., i^ref) with = 0, but 
we must determine the relative angles (p with respect to the reference Hamiltonian for all other 
control settings. 



4.1. Rotation frequency and declination of rotation axis 

To relate the Hamiltonian parameters to observable dynamics, it is instructive to visualize qubit 
states and their evolution on the Bloch sphere. If we define s = {sx, Sy, s^)^ with Sk = Tr(/)(Tfc) 
for k E {x, y, z} then it is easy to check that there is a one-to-one correspondence between 
density operators p and points s inside the closed unit Ball in IR^. Furthermore, the evolution 



of p{t) under the (constant) Hamiltonian (12) corresponds to a rotation of s{t) about the (unit) 
axis d = (sin cose/), sin sin (/>, cos ^)"^ with angular velocity to as illustrated in Fig. |4| The 
figure also shows that we can determine the angle 9 and angular frequency ||d|| = uj from the 
projection of s(t) under the rotation about the axis d onto the z-axis. As the figure shows, given 
z{t), we can in principle extract lo and 9 from the first minimum (to,-Zo) of -2^(0) although in 
practice it is usually preferable to use Fourier analysis or harmonic inversion techniques. z(t), 
being the expectation value {A{t)) = Tv[azp{t)] of the observable A = az m the state p{t), can 
be obtained experimentally as follows (|12p : 



(i) Initialize: Measure and record outcome Aq = ±1 ^ system in state AqSq = Aa(0, 0, 1)'^. 

(ii) Evolve: Let the system evolve for time t under the control settings f =^ system now in 
state Aas(t) = Xa{x{t),y{t), z{t))'^ with 

z{t) = cos^ 9 + sin^ 9 cos(u;t) . (13) 

(iii) Measure: Repeat measurement and record outcome A^ = ±1 =^ system now in state A^sq. 

(iv) Repeat: steps (i)-(iii) A^e times. 

Let Na^b be the number of times the the initial measurement produced outcome Aa and the final 
measurement produced outcome A;,. Then there are four possible combinations of outcomes 
-^0,0) -^0,1) -^1,0 and A'^1,1, which must add to the total number of experiments Ne- The number 
of experiments that started in the state |0) [or sq = (0,0,1)"^] is = Nq^ + A''o,i, and the 
fraction of experiments for which the second measurement yields Aq conditioned on the system 
starting in the state |0) is Nofi/N^. Thus the ensemble average of ^ = o"^ at time t, assuming 
we started in the state |0), is 

, Nofl , , iVo,i AoiVo.o + AiiVo.i 

^°lvr + ^^lvr- iVo,o + iVo,i ' ^''^ 

and for Ng ^ 00 this relative frequency should approach the true expectation value of z{t) 
assuming z(0) = 1. Thus, we can in principle determine z{t) to arbitrary accuracy by choosing 
Ne large enough. By repeating the experiments for different evolution times tk, e.g., using a 
stroboscopic mapping approach (see Fig. [5]), we can determine z{t) as a function of t, from which 
we can extract 00 = ujf and 9 = 9f as discussed. For instance, noting that z{t) = 2pQ{t) — 1, Fig.js] 
shows that the 1st order Fourier peak F{lij) of z{t) in this example is F{uj) « 2 x 0.2468 = 0.4935 



for (jj ~ 2 and thus comparison with Eq. (13) shows that sin^ 9 = 2^[F{uj)], or equivalently, 
9 ~ 1.4566. These values are reasonably good estimates for the real values uj = 2.0086 and 
9 = —1.4780, except for the sign of 9, or the orientation of the rotation axis, which we cannot 
determine from the given data, as these rotations have identical projections onto the z-axis. 



This is reflected in Eq. (13) by the fact that both coefficients cos 9 and sin 9 contain squares. 
In principle, we can determine the parameters uj and 9 to arbitrary accuracy using this approach 
based on regular sampling and Fourier analysis (|13p . although the total number of experiments 
can become prohibitively large. An alternative that merits further investigation is the use of 
adaptive sampling techniques to reduce the total number of experiments necessary. 



4-2. Relative angles between rotation axes 

Having determined the parameters ujf and 9f for a particular control setting f and chosen a 
suitable reference Hamiltonian i/j-ef with (pj-ef = 0, to complete the characterization of Hf 7^ i/ref , 
we must determine the horizontal angle (pf. The reference Hamiltonian H^e{ must not commute 
with the measured observable A, or equivalently the angle ^ref between the rotation axis dref 
and the z-axis must be nonzero. Ideally dref should be as close to orthogonal to the z-axis as 
possible. Assuming ^ref S (7r/4, 7r/2], and ujf and 9f are known, we can determine the angles (pf 
by performing the following experiments (12): 

(i) Initialize: Measure and record outcome Aa = ±1 ^ system in state AqSq = Aa(0, 0, 1)'^. 

(ii) Prepare: Rotate around reference axis d^ef by angle uq = arccos(j^), x = — cos(20ref) 
^ system in new state AaSi = Aa(cos (3, sin /?, 0)"^ with j3 = arctan(— -v/22;) / cos ^ref • 




Figure 5: Stroboscopic mapping of coherent oscillations with projective measurements. After 
initialization through measurement, the system is allowed to evolve for fixed times under the 
influence of the controls before a second measurement is taken. The experiments are repeated 
to determine the expectation value of the observable at each time. 



(iii) Evolve: Let the system evolve for time t under the control settings f 
state Aas(t) = Xa{x{t),y{t), z{t))'^ with 

z{t) = c[l — cos(cjft)] + dsin(a;ft) 
c = sin 6f cos Of cos((/)f — /?) 
d = sin 6{ sin(0f — /?) . 
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(iv) Measure: Repeat measurement and record outcome = ±1 system in new state A^sq. 

(v) Repeat: steps (i)-(iv) A'^e times. 

As before, if Na^ is the number of times the initial measurement produced outcome Aa and the 
final measurement produced outcome Af, then we have 



z{t) 



XpNofl + Ai Ao,i 



for sufficiently large Ng, thus allowing us to determine z{t) experimentally. Since 9f, ujf and /3 are 



known from the 9, u> characterization step, Eq. ( 15 ) allows us to determine (f>f via the coefficients 



c and d, which can be determined in principle either through curve fitting or by taking the 
Fourier transform of z(t) and noting that c and d correspond to the real and imaginary part of 
the 0th and 1st order Fourier peaks. 

As a specific example, assume we have two Hamiltonians Hi and H2 and we have already 
determined uJi = 1, 9i = j, and a;2 = 1.2, ^2 = f • We choose Hi as reference Hamiltonian and 
note that a rotation about Hi hy a = tt maps the state sq = (0,0,1)^ to si = (1,0,0)^. We 
then let state si evolve for different length of time t under Hamiltonian H2, and determine the 



projection z(t). This will produce a measured z-trace as shown in Fig. pi Eq. (15) shows that 
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Figure 6: Determination of relative horizontal angle between two Hamiltonians. Measured 
projection z{t) of si(t) onto the z-axis and corresponding Fourier spectrum. Blue/red 
bars indicate the cosine/sine coefficients, i.e., real/imaginary parts of the Fourier transform, 
respectively. 



in theory the height of the real part of the 0th and 1st order Fourier peaks should be the same 
except for the sign, but the graph clearly shows that this is not the case for our noisy (simulated) 
experimental data. In fact, the 1st order peaks around itwf are significantly broadened, while 
the 0-peak is very narrow. This suggests that we are likely to obtain the best estimate for <j) 
by setting c = FT[z](0) and (pf = arccos[c/(sin 6'f cos6'f)] + /3, where /3 = in our case. Indeed, 
this yields an estimate of (pf = 0.25187r, which is very close to actual value tt/4 used in the 
simulation. Note that using the 1st order Fourier peaks gives far less accurate results due to 
the significant broadening of the peaks. This can be reduced using phase matching conditions 
and other tricks (13), but in this case simply using only the 0-peak is sufficient to get a good 
estimate. 



4-3. Functional dependence on control settings 

In some settings it is sufficient to characterize the Hamiltonians Hf for a discrete, finite set 
of controls {f^'^)}. This would be the case for bang-bang-type control schemes that require 
only switching between a finite set of control Hamiltonians. In other cases, in particular 
when continuously varying fields are to be applied, we would also like to know the functional 
dependence of the Hamiltonian on the control fields, e.g., whether we can assume a control-linear 
model as in (|4]), or if there are non-linear or crosstalk effects, etc. Although characterizing the 
control dependence of the Hamiltonian is a non-trivial problem, given sufficiently many Hf data 
points, it is possible to test the quality of fit of a particular model using statistical means. For 
instance, suppose we have characterized the Hamiltonians Hf for f G {f^'^^}. Let be the 

corresponding vectors in IR as defined above, and let fm be the mth component of the control 
vector f^'^). We can find the best-fitting linear model by finding vectors G IR^, m = 0, . . . , M, 
where M is the number of independent control variables, that minimize the norm of the residuals 



0.7 















Li -Li 


E J 


Lrj>rVJlll.rLiI.JLJI.£i..'LJlJJILJLJLIkB 






w ^ 


V 






















, , evolution time t, 


N =100, dt=0.01, f=0.6 

, «: ■ 1 



10 20 30 40 50 00 70 80 90 1 0O 




control inputf farb. unrts) 



Figure 7: Determination of functional dependence of Hf on control input f . The bottom graph 

(k) (k) 

shows the estimated values for dx and dl (circles and triangles, respectively) as well as the 
actual functional dependence in the model (solid lines). The top graph shows a typical z{t) 
trace obtained by averaging over = 100 simulated experiments for each time = O.Olfc. The 
effects of finite sampling and projection noise are clearly visible for this trace. Nonetheless, the 
predicted functional dependence of Hf on the control input / is in good agreement with the 
actual dependence in the model. 



^ = yEfc llrfclli where 

M 

rfc = do + ^ fi^dm - dk. (16) 



m=l 



The norm of the residuals A then gives an indication of the quality of fit of a linear dependence 
model, and allows us to compare the quality of fit for different models. 

As an example, consider a Hamiltonian Hf depending on a one-dimensional control input 
f = /. To characterize the functional dependence of Hj on /, we choose a finite set of control 
settings fk = 0.2(k — 1) for k = 1,...,10. We characterize first the rotation frequencies uJk 
and declination angles 6^ for each f^, and then the relative angles (pk as discussed above. From 



these, we calculate the components dl^\ d^"* and d^z^ of the Hamiltonian according to Eq. (|12|), 



and plot them versus / as shown in Fig. [7] Specifically, in our case we find (f)^ ~ and hence 
dy = 0, dz appears to be constant and dx exhibits a quadratic dependence on the control /. 
A least square fit of the data using a linear and quadratic fitting function, respectively, yields 
dz{f) = 1.002 - 0.003/ and dx{t) = 0.006/ + 0.996/^, which is in very good agreement with 
the actual functional dependence in the model dz[f) = 1 and dx{f) = /^ despite the fact that 
the simulated data appears very noisy and the 9 and oj estimates did not employ tuning such as 
phase matching to further improve the accuracy. 



4-4- Characterization in the presence of decoherence 

The protocols for Hamiltonian identification above are sufficient if dissipative effects are small on 
the timescales for coherent control, i.e., if the system and environment coupling is much weaker 



than the system-control interaction. If the couphng between the system and its environment is 
stronger then dissipative effects must be directly incorporated into the characterization protocols. 
Assuming Markovian decoherence of the form ([2|, it is easy to show that for = 2, we 



that defines a preferred decoherence basis. The operator (Tz is usually associated with pure 
decoherence (pure dephasing) while the raising and lowering operators (t± are associated with 
population relaxation from |1) to |0) and vice versa. Very often it is furthermore assumed that 
U = I is the identity, i.e., that the preferred basis for relaxation is the same as the measurement 
basis. Although this need not be the case, it is often a good approximation. The main effect 
of relaxation processes is to dampen the observed coherent oscillations z(t) that result from the 
Hamiltonian dynamics, and try to force the system 'asymptotically' into a 'steady' state. This 
effect manifests itself in Lorentzian broadening of the peaks in the oscillation spectrum, which 
can be used for characterization purposes. 

For example, Fig.[8]shows significantly damped oscillations as a result of relaxation. The fact 
that the oscillations decay to Zoo = suggests either pure dephasing, or a symmetric relaxation 
process with T^ = r_, i.e., equal probability of relaxation from state |0) to |1) and vice versa, 
or a combination of both. Since symmetric relaxation and pure dephasing have very similar 
signatures, it can be difficult to distinguish these processes from the oscillation trace alone. 
However, if there exists a control f such that Hf oc a^, which corresponds to a rotation about 
the z (or measurement) axis, then we can in principle distinguish the two processes by initializing 
the system repeatedly in either state |0) or |1) and letting the system evolve under Hf for various 
periods of time before measuring again. If the decoherence process is pure dephasing then the 
system will remain in the initial state, while we expect to observe random jumps (with equal 
probability) for symmetric relaxation. Assuming this preliminary characterization step suggests 
that the damping in the figure is due to pure dephasing, we can then extract information about 
the rotation frequency lvq and dephasing rate T by fitting a Lorentzian envelope function 



to the first-order peak in the Fourier spectrum as shown in the figure. In our example, using 
a simple least-squares minimization yields estimated values for both ujq and F that are close to 
the actual values used in the simulated experiment. This approach can be generalized for more 
complicated relaxation processes (jl4|] . Alternatively, if we can characterize the Hamiltonian 
dynamics sufficiently to be able to initialize the system in non-measurement basis states and 
simulate readout in different bases, then generalized Lindblad operators can be estimated using 
repeated process tomography (ITSll . 

5. Discussion of generalizations and conclusions 

We have described some basic tools and techniques for characterizing the extent to which 
the dynamics of a certain system is confined to a low-dimensional subspace of a potentially 
large Hilbert space, as well as for identifying both Hamiltonian and decoherence parameters 
experimentally, with emphasis on protocols that are realistic even in the presence of limited 
direct readout and control capabilities. 

Although the explicit schemes presented focused on qubit-like systems, the same basic ideas 
can be applied to higher-dimensional systems. However, the number of Hamiltonian and 
decoherence parameters to be determined for higher- dimensional system makes finding explicit 
protocols for complete characterization quite challenging. The qubit characterization protocols 
show that identifying the relative angles between several control Hamiltonians requires a two- 
stage process and more complicated two-step experiments (initialization in non-measurement 
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Figure 8: Characterization in the presence of decoherence. The observed damped oscihations 
due to pure dephasing (top) result in Lorentzian broadening of the Fourier peaks, which can be 
exploited to obtain (good) estimates for the rotation frequency and decoherence rate F. 



basis state followed by controlled evolution) than identifying the rotation frequency and 
declination angle for a single Hamiltonian. As the dimension of the system increases more 
parameters are necessary to fully characterize the Hamiltonians, leading to more complex multi- 
stage characterization protocols. 

Another problem is that the number of experiments required to implement such schemes 
may become prohibitively large for parameter estimation based on spectroscopic analysis. One 
of the most promising strategies to avoid this problem is the development of efficient adaptive 
characterization schemes — as opposed to protocols based on regular sampling — to minimize the 
number of measurements / experiments required. For systems composed of separate smaller units 
such as multi-qubit systems, boot-strapping approaches to characterization may also provide a 
fruitful alternative. For example, assuming we have separately characterized individual qubits 
using qubit identification protocols and therefore have full local control, we can reduce the 
problem of two-qubit identification in principle to identifying the interaction Hamiltonian. For 
two-qubit systems with local control and a fully non-local interaction Hamiltonian this can be 
achieved using entanglement mapping or concurrence spectroscopy (\TU\ llZp . 
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